#include <iostream>
#include<sys/time.h>
#include<omp.h>
#define N 4800
#define num_thread 7
using namespace std;
int n=100;
float m[N][N],init_m[N][N];
void init(int n){
    for(int i=0;i<n;i++){
        for(int j=0;j<n;j++)
            init_m[i][j]=0;
        init_m[i][i]=1.0;
        for(int j=i+1;j<n;j++)
            init_m[i][j]=rand()%1000+1;
    }
    for(int k=0;k<n;k++)
        for(int i=k+1;i<n;i++)
            for(int j=0;j<n;j++)
                init_m[i][j]=int((init_m[i][j]+init_m[k][j]))%1000+1.0;
}
void reset(int n){
    for(int i=0;i<n;i++)
        for(int j=0;j<n;j++)
            m[i][j]=init_m[i][j];
}
//平凡算法
void serial(int n){ //串行
    struct timeval h, t;
    float time = 0.0;
    gettimeofday(&h, NULL);
      for(int k=0;k<n;k++){
          for(int j=k+1;j<n;j++)
              m[k][j]/=m[k][k];
          m[k][k]=1.0;
          for(int i=k+1;i<n;i++){
              for(int j=k+1;j<n;j++)
                  m[i][j]-=m[i][k]*m[k][j];
              m[i][k]=0;
          }
      }
    gettimeofday(&t, NULL);
    time += 1000*(t.tv_sec - h.tv_sec) + 0.001*(t.tv_usec - h.tv_usec);
    cout<<"serial: "<<time<<" ms"<<endl;
    return;
}
//openMP
void om() {
    struct timeval h, t;
    float time = 0.0;
        gettimeofday(&h, NULL);
        int i,j,k;
        float tmp;
        #pragma omp parallel num_threads(3),private(i, j, k, tmp)
        for(k = 0; k < n; k++) {
            #pragma omp master
            {
                tmp = m[k][k];
                for(j = k + 1; j < n; j++) {
                    m[k][j] = m[k][j]/tmp;
                }
                m[k][k] = 1.0;
            }
            #pragma omp for schedule(guided, 5)
            for(i = k + 1; i < n; i++) {
                tmp = m[i][k];
                for(j = k + 1; j < n; j++) {
                    m[i][j] = m[i][j] - tmp*m[k][j];
                }
                m[i][k] = 0;
            }
        }
        gettimeofday(&t, NULL);
        time += 1000*(t.tv_sec - h.tv_sec) + 0.001*(t.tv_usec - h.tv_usec);
    cout<<"om: "<<time<<" ms"<<endl;
}
void om1() {
    struct timeval h, t;
    float time = 0.0;
        gettimeofday(&h, NULL);
        int i,j,k;
        float tmp;
        #pragma omp parallel num_threads(7),private(i, j, k, tmp)
        for(k = 0; k < n; k++) {
            #pragma omp master
            {
                tmp = m[k][k];
                for(j = k + 1; j < n; j++) {
                    m[k][j] = m[k][j]/tmp;
                }
                m[k][k] = 1.0;
            }
            #pragma omp for schedule(guided, 5)
            for(i = k + 1; i < n; i++) {
                float tmp = m[i][k];
                for(j = k + 1; j < n; j++) {
                    m[i][j] = m[i][j] - tmp*m[k][j];
                }
                m[i][k] = 0;
            }
        }
        gettimeofday(&t, NULL);
        time += 1000*(t.tv_sec - h.tv_sec) + 0.001*(t.tv_usec - h.tv_usec);
    cout<<"om1: "<<time<<" ms"<<endl;
}
void om2() {
    struct timeval h, t;
    float time = 0.0;
        gettimeofday(&h, NULL);
        int i,j,k;
        float tmp;
        #pragma omp parallel num_threads(7),private(i, j, k, tmp)
        for(k = 0; k < n; k++) {
            #pragma omp master
            {
                tmp = m[k][k];
                for(j = k + 1; j < n; j++) {
                    m[k][j] = m[k][j]/tmp;
                }
                m[k][k] = 1.0;
            }
            #pragma omp for schedule(guided)
            for(j = k + 1; j < n; j++) {
                float tmp = m[k][j];
                for(int i = k + 1; i < n; i++) {
                    m[i][j] = m[i][j] - tmp*m[i][k];
                }
            }
            for(i = k + 1; i < n; i++) {
                m[i][k] = 0;
            }
        }
        gettimeofday(&t, NULL);
        time += 1000*(t.tv_sec - h.tv_sec) + 0.001*(t.tv_usec - h.tv_usec);
    cout<<"om2: "<<time<<" ms"<<endl;
}

int main()
{
cout<<"n="<<n<<endl;
    init(n);
    reset(n);
    serial(n);
    reset(n);
    om();
    reset(n);
    om1();
    reset(n);
    om2();
    n=200;
    cout<<"n="<<n<<endl;
    reset(n);
    serial(n);
    reset(n);
    om();
    reset(n);
    om1();
    reset(n);
    om2();
    n=400;
    cout<<"n="<<n<<endl;
    reset(n);
    serial(n);
    reset(n);
    om();
    reset(n);
    om1();
    reset(n);
    om2();
    n=800;
    cout<<"n="<<n<<endl;
    reset(n);
    serial(n);
    reset(n);
    om();
    reset(n);
    om1();
    reset(n);
    om2();
    n=1600;
    cout<<"n="<<n<<endl;
    reset(n);
    serial(n);
    reset(n);
    om();
    reset(n);
    om1();
    reset(n);
    om2();
    n=3200;
    cout<<"n="<<n<<endl;
    reset(n);
    serial(n);
    reset(n);
    om();
    reset(n);
    om1();
    reset(n);
    om2();
    n=4800;
    cout<<"n="<<n<<endl;
    reset(n);
    serial(n);
    reset(n);
    om();
    reset(n);
    om1();
    reset(n);
    om2();
    return 0;
}
